Theoretical Analysis of the "Double-q" Magnetic Structure of CeAl2 
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A model involving competing short-range isotropic Heisenberg interactions is developed to explain 
the "double-q" magnetic structure of CeAl2. For suitably chosen interactions terms in the Landau 
expansion quadratic in the order parameters explain the condensation of incommensurate order at 
wavevectors in the star of (1/2 — 5, 1/2 -|- 5, l/2){2n/a), where a is the cubic lattice constant. We 
show that the fourth order terms in the Landau expansion lead to the formation of the so-called 
"double-q" magnetic structure in which long-range order develops simultaneously at two symmetry- 
related wavevectors, in striking agreement with the magnetic structure determinations. Based on 
the value of the ordering temperature and of the Curie- Weiss G of the susceptibility, we estimate 
that the nearest neighbor interaction, Kq, is ferromagnetic with Ko/k = — 11 ± IK and the next- 
nearest neighbor interaction J is antiferromagnetic with J/k = 6 ± 2K. We also briefly comment on 
the analogous phenomena seen in the similar system, TmS. 

PACS numbers: 75.25.-(-z, 75.10.Jm, 75.10.Dg 
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I. INTRODUCTION 

CeAl2 (CEAL) is a metallic system whose magnetic 
structure has been the object of some controversy for 
several yearsi^ii^i^ Initial studiesii^ indicated the ex- 
istence of incommensurate long-range magnetic order on 
the Ce ions with a single wavevector in the star of q, 
where 

q = (2n/a){l/2- 6,1/2 + 6,1/2) = {2Tr/a)q, (1) 

with 6 = O.ll.ii^ Later— it was proposed that this struc- 
ture involved the simultaneous condensation of three 
wavevectors in the star of q, but this suggestion of a 
"triple q" structure was refuted in Rcf. 4. More recent 
worli^ showed that the structure was in fact a "double 
q" one in which exactly two wavevectors in the star of 
q were simultaneously condensed. In addition, contin- 
ued interest in CEAL is due to its Kondo-like behavior. 
Initial indications of this came from the observation of a 
minimum in the resistivity at about 15K, which was at- 
tributed to spin compensation.^ A single impurity model, 
with a Ce"^+ ion in a cubic crystal field interacting with 
the conduction band, was able to account for most of 
the electrical properties.S Moreover, when, in neutron ex- 
periments, no third order magnetic satellite appeared at 
low temperature, the Kondo effect was invoked to ex- 
plain why the moment of a Kramers ion did not sat- 
urate in the zero-temperature limit. ^ This objection is 
partially removed by the double-q structure.^ Moreover, 
an analysis^ of multi-q states claims that the double-q 
structure can not be explained if CEAL is regarded as 
an itinerant-electron magnet. 

In this paper we proceed under the assumption that 
although Kondo effects may be present due to the cou- 
pling of the Ce 4f electron to the conduction band, the 
magnetic structure can be understood in terms of inter- 
action between localized moments on the Ce ions. Since 
the lattice structure is fee, it is apparent that antifer- 



romagnetic interactions between shells of near neighbors 
could compete and might then explain the incommensu- 
rability. However, no concrete calculations of this type 
have yet appeared. It is also interesting that this system 
does not follow the simplest scenario for incommensu- 
rate magnetSfiSi namely, as the temperature is lowered, 
a phase transition occurs in which a modulated phase 
appears with spins confined to an easy axis, and then, 
at a lower temperature a second phase transition oc- 
curs in which transverse order develops, so as to par- 
tially satisfy the fixed length spin constraint expected 
to progressively dominate as the temperature is lowered. 
Instead, in CEAL, there is no second phase transition, 
and in the ordered phase one has the simultaneous con- 
densation of long-range order at two symmetry-related 
wavevectors.^'- There are two aspects of this behavior 
that have not yet been explained. 1) The incommen- 
surate wavevector lies close to, but not exactly along 
the high symmetry (1, 1, 1) direction and 2) although so- 
called "triple-q" systems are well known yii^i^iiiiiii^ in 
which the incommensurate ordered state consists of the 
simultaneous superposition of three wavevectors, it is un- 
usual, in a cubic system, to have a "double-q" stateiSiii 
consisting of the simultaneous superposition of exactly 
two wavevectors. 

The aim of this paper is to develop a model which can 
explain the above two puzzling features. We first ad- 
dress the determination of the incommensurate wavevec- 
tor. Some time ago, Yamamoto and Nagamiya^^ (YN) 
studied the ground state of a simple fee antiferromag- 
net with isotropic nearest-neighbor (nn) and next-nearest 
neighbor (nnn) Heisenberg interactions and found a rich 
phase diagram in terms of these interactions whose cou- 
pling constants we will denote here as J and M, respec- 
tively. We perform an equivalent calculation for a re- 
lated model appropriate to CEAL based on an analysis 
of the terms in the Landau expansion of the free energy 
in the paramagnetic phase. By studying the instability 



of this quadratic form which occurs as the temperature 
is lowered, one can predict the magnetic structure of the 
ordered phase. In particular, one can thereby determine 
the wavevector at which this instability first occurs. This 
phenomenon is referred to as "wavevector selection." As 
Nagamiya's reviewsiS indicates, correct wavevector se- 
lection in CEAL must require a model which involves 
competition between nn and further neighbor interac- 
tions. For the fee structure of CEAL the most con- 
venient model which almost explains wavevector selec- 
tion involves nn, nnn, and fourth-neighbor interactions. 
Based on our insight developed from this model, we sug- 
gest the how more general interactions can completely ex- 
plain wavevector selection. Although we invoke more dis- 
tant than nn interactions, the magnitudes of the further 
neighbor couplings needed to explain the nonsymmetric 
wavevector of CEAL decrease with increasing separation 
and are reasonable, especially in view of the possibility of 
Ruderman-Kittel-Kasuya-Yosida (RKKY) ^^ interactions 
in this metallic system. Because our main interest lies 
in explaining wavevector selection, we have completely 
ignored anisotropy, whose major effect is to break ro- 
tational invariance and select spin orientations. Coinci- 
dentally we note several regions in parameter space for 
these models in which one has a multiphase point (at 
which wavevector selection is incomplete). This phe- 
nomenon is perhaps most celebrated in the Kagome^S 
and pyrochloreSi systems. Based on these results we also 
point out that there are likewise regions of parameter 
space that could explain wavevector selection22i2i24 in 
the similar Kondo-like system TmS. 

The second stage of our calculation for CEAL involves 
an analysis of the fourth order terms in the Landau ex- 
pansion, because it is these terms which dictate whether 
only one or more than one wavevector in the star of q 
is simultaneously condensed to form the ordered phase. 
For this analysis there are two plausible ways to proceed. 
An oft-used approacbii is to determine the most gen- 
eral fourth order term allowed by symmetry and then see 
whether some choice of allowed parameters can explain 
a "double q" state. The virtue of this method is that it 
corresponds to the use of fluctuation-renormalized mean- 
field theory. A drawback, however, is that it is hard to 
know whether the allowed parameters are appropriate for 
the actual system. Here we adopt a contrary procedure 
in which only the "bare" (unrenormalized) fourth order 
terms are considered. Obviously, these terms do have 
the correct symmetry, and although they might not be 
the most general possible terms, they do ensure that the 
values of the parameters are plausible. 

The organization of this paper follows the above plan. 
In Sec. II we extend the analysis of YN to fee magnets 
with three shells of isotropic exchange interactions, but 
even this model only partially explains the wavevector 
selection seen in CEAL. In Sec. Ill we invoke more dis- 
tant interactions, whose existence is attributed to either 
RKKY interactions,^^ or indirect interactions via excited 
crystal field states, as discussed in Appendix C. Thereby 



we explain wavevector selection in CEAL and also in the 
similar system TmS. Here we also use the observed or- 
dering temperature and data for the zero wavevector sus- 
ceptibility to estimate values of the dominant exchange 
interactions. In Sec. IV we analyze the fourth order 
terms in the Landau expansion and show that they nat- 
urally lead to the "double-q" state observeci§i^ in CEAL. 
Our results are briefly summarized in Sec. V. 



II. WAVEVECTOR SELECTION FOR A "3- J" 
MODEL 

In isotropic Heisenberg models of magnetic systems 
with only nearest neighbor (nn) interactions on, say, a 
simple cubic lattice, the magnetic structure of the or- 
dered phase is trivially constructed if the sign of the inter- 
action is known. In more complicated models it may hap- 
pen that next-nearest neighbor (nnn) interactions com- 
pete with the nn interactions, in which case the mag- 
netic structure may be an incommensurate onei^ In this 
case, the quadratic terms in Landau free energy (which 
we study below) will be such that, as the temperature is 
lowered, the paramagnetic phase develops an instability, 
relative to the development of long-range magnetic or- 
der, at a wavevector q (or more properly, at the star of 
q) . For CEAL our aim is to study this "wavevector selec- 
tion," and explain how a model of exchange interactions 
can lead to the observed ordering wavevectors. 

For this purpose, this section is devoted to an analysis 
of the quadratic terms in the free energy which determine 
wavevector selection. We first note that the space group 
of CEAL is Fd3m (space group #227 in Ref. HI) which 
is an fee system with two Ce atoms per fee unit cell at 
locations 



Ti- (0,0,0), X2 = (l,l,l)(a/4) 



(2) 



This means that each Ce ion has a tetrahedron of nn's 
and we will treat further neighbor interactions as in an 
fee Bravais lattice. Note that the two sites at Ti and T2 
are related by inversion symmetry relative to the point 
(1, 1, l)(a/8). We introduce the following simple model 
of exchange interactions, 

n = Y. 4°,U^'^')Sop(R + x„).Sop(R' + x„(^) 

R.,n;R,',n' 

where Sop(R.-t- t„) is the spin operator at R + t„. Here 
we treat the model having three shell of interactions, so 
that the only nonzero J's are 

Jl2 (R' R ) = -^21 R J R) 

= X(°) if |R + Ti-R'-T2| = a\/3/4 

ji?R,R') = 42^R,R') = J^•'^ if |R-R'|=a/V2, 

j}?R,R') = 42'R,R') = M^°^ : if |R- R'l = a . (4) 

In other words we have exchange couplings, K'-^\ J(°\ 
and Af(°) between nn's, nnn's, and a shell of fourth- 
nearest neighbors (fnn's), respectively, and these are 
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FIG. 1: (Color online) Exchange interactions for CEAL, 
showing only the magnetic Ce atoms in three (001) planes. 
The structure can be thought of as consisting of two inter- 
penetrating fee sublattices. The Ce t\ atoms are circles (red 
online) and the Ce T2 atoms are squares (blue online). Filled 
squares represent the four Ce T2 atoms which are nearest 
neighbors to the Ce t\ atom at the origin, indicated by a 
filled circle. The exchange constants between nn's, nnn's, 
third neighbors, and fourth neighbors are K, J, L, and A/, 
respectively. 



shown in Fig. ^ Equation (^J implies that positive ex- 
change constants are antiferromagnetic. Our j'"' and 
jVf (0) correspond to YN's — Ji and —J2, respectively and 
their simple fee structure did not have a K^^^^ interaction. 
As will become apparent below, we include a fnn inter- 
action rather than a third neighbor (tnn) interaction in 
the interest of algebraic simplicity. 

The 4f electron of the Ce ion has quantum numbers 
L^3, S ^ 1/2, and J = |L -h S| = L - 5 = 5/2, so that 
Sop — {g,j — 1) J, where gj = 6/7 is the Lande g factori^ 
The crystal field then splits the six states of the J = 5/2 
manifold into a ground doublet and an excited quartet 
state at an excitation energy of about lOOK in tempera- 
ture units i2iSL22iSS This ground doublet can be described 
by an effective spin operator Sos of magnitude 1/2 and 
within the doublet J = (5/3)Scff when admixtures from 
the quartet state are neglected. In that case 



Sop = (5/3) (gj - l)Soff = .goScff 



(5) 



When admixtures caused by the actual exchange field 
and also an applied field of 45 kOe were calculated by 
Barbara et al.,'^'^ the moment was found to be somewhat 
larger than that zero net field, but for zero applied field 
we neglect this effect. Then we write the Hamiltonian in 
terms of effective spins 1/2 as 



n = 



/ ^ J„^„'(R, R')Soff(R + T„) • Soff(R' 



R.,n;R,',n 



where 



Jnn' (R, R ) — g^Jrin' (^' ^') ' 



2 r(o) 



(7) 



and we have the interactions K , J, and M analogous to 
those in Eq. (gj). 

We now develop the Landau expansion for the free en- 
ergy. The approach we follow is to write the trial free 
energy as 



F = Tr[p7^-|-fcrplnp] 



(8) 



where p is the trial density matrix which is Hermitian and 
has unit trace. The actual free energy is the minimum of 
F with respect to the choice of p. Mean field theory is 
obtained by restricting p to be the product of single-spin 
density matrices, so that 



P = Wpi^-^n) 



(9) 



R.ri 



where p(R, t„) is the density matrix for the Ce spin at 
R -H T„. We write 



p(R + -r„ 



1 



[l + a(R + T„).S(R + T„)] ,(10) 



where from now on S(R + t„) denotes the effective spin 
1/2 operator for the site in question and we identify the 
vector trial parameter a by relating it to the thermal 
expectation value of the spin as 

(S(R-Ht„)) = Tr[p(R + T„)S(R + T„)] 

= a(R + r„)/4, (11) 

so that 

p(R + x„) = i[l + 4(S(R,T„))S(R,x„)] . (12) 
Then, one finds that 

F = \ Y. >/«,n'(R,R')(S(R + x„)).(S(R' + T„0) 
-TS , (13) 



where 



TS = A:T^TrUh+4(S(R,r„)>S(R,T„) 

R,n I \ 



xln -[l + 4(S(R,T„))S(R,x„)] ,(14) 



which we evaluate as 



R,n p—1 



[4(5(R + T„))-(5(R + r„))]P 



2p{2p-l) 



(15) 



In this section we consider the term quadratic in the spin 
variable and in the next section we consider the quartic 
term in this expansion. (Higher order terms are not nec- 
essary for our analysis.) 

We introduce as order parameters, the Fourier coeffi- 
cients defined for n = 1, 2 by 

(S(R + T„)) = S„(q)e^'i-^ + 5„(q)*e-'i-^. (16) 

Note that the phase factor is determined by the origin 
of the unit cell and not by the actual location of the 
spin site. Any two wavevectors which differ by a linear 



combination of reciprocal lattice basis vectors, G„ are 
equivalent, where 



Gi = 2^(l,l,-l)/a 
G2 = 2^(l,-l,l)/a 
G3 - 2^(-l,l,l)/a 



(17) 



Now we write the contribution to the free energy which 
depends on the order parameter for some wavevector q. 
In terms of this order parameter, the mean field free en- 
ergy at quadratic order, F2, can be written as 



^2 = lj2^X-X,n'Sn{cirSn>{q) , 



where^ 

x-'(q) = 

where 



4fcT + Jn(q) Ji2(q) 

J2i(q) 4fcT + J22(q) 



(18) 



(19) 



<^ii(q) = 2M{cosqa;a-\- cosqya -{- cosq^a) 

= ^22 (q) 

Ji2(q) = ii'[l + e-'(9-+«'')''/2 ^ e-'('?x+g.)a/2 

^^-^{1y+q.)a/2^ ^ J2l(q)* , (20) 

where Cq. = cos(gcta/2). Omitting the factor 4fcr, the 
minimum eigenvalue of the x~^ matrix, which selects 
the wavevector, is 



A(q) = Jll(q)-|Jl2(q)| 

= 2Af (2c^ + 2c2 + 2c2 - 3) 

+AJ{R^ -1)-2\K\R , 



where 



R 



[1 -|- CxCy -j- CxCz \ CyCz\ 



1/2 



(21) 



(22) 



(By the square root, we always mean the positive square 
root.) Note that changing the signs of all the c's corre- 
sponds to adding a reciprocal lattice vector to q and does 
not change A(q). So solutions which differ by changing 
the signs of all the c's are equivalent to one another. Al- 
though the minimum value of the free energy does not 
depend on the sign of K, the ratio of spin amplitudes 
within the unit cell does depend on this sign. To discuss 
the sign of K it is convenient to set q = (1, 1, l)(7r/a) 
(which is nearly the wavevector of interest). Then if K 
is negative (ferromagnetic), Ji2(q) is negative and the 
minimal spin eigenvector is (1,1), which indicates that 
the spins at xi and T2 are parallel, as is illustrated in 
Fig. 13 whereas if K is positive, they are antiparallel. In 
the former (latter) case, the other three spins of the nn 
tetrahedron are antiparallel (parallel) to the spin at the 
origin. Thus the sign of K is easily related to whether 
the majority of the nn's are parallel in which case K is 




-a/4 
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FIG. 2: (Color online) Nearest neighbor interactions K (in- 
dicated by dashed lines). Here we show planes perpendicular 
to (111), with Ce(l) sites represented by circles (red online) 
and Ce(2) sites represented by squares (blue online). Below 
each plane the value oi x + y + z for that plane is given. If the 
wavevector is assumed to be {n, n, n)/a, then all spins within 
a given plane are parallel to one another and if K is negative 
(ferromagnetic) , then the spin directions of each plane are as 
indicated by the arrows. Had K been of opposite sign, then 
the spin directions of the Ce(2) would be reversed and an 
inequivalent magnetic structure would be realized. 



negative. Otherwise K is positive. The structure deter- 
minations indicate that the correct choice is that K is 
negative (ferromagnetic). Henceforth K will be used to 
denote \K\. 

As mentioned in the introduction, this system for K — 
has been comprehensively analyzed by YN. However, 
they seem to have overlooked an amusing limit for K = 0. 
Namely, if J = 2Af we have 

A(q) = -3J + 2J{c^ + Cy+c^)^ . (23) 

For J < this is minimal for Cx = Cy = Cz = H. For 
J > this is minimized over the entire two dimensional 



manifold for which Cx -t- Cy 



0. What this means 



is that for this special case, there is no wavevector selec- 
tion. Such a multiphase point has been found in several 
modelsi^SiSLS^ As we shall see, this multiphase behav- 
ior is modified to encompass a one dimensional manifold 
when K is small and M < J/ 2. 



A. J<0 

When J < 0, then the second and third terms of Eq. 
(|21|1 do not compete with one another: for a fixed value 
oi c^ + c^+cl, A(q) is minimized by maximizing R, which 



implies that c^ 



so that 



A(q) = ~6M + 12( J -I- M)c^ - 2Ky/l + 3c2 . (24) 



The extrema must be either for (? 
differentiation) 



3c^ 



-1 



K 



_4(J + M) 

so that, for this to apply, we must satisfy 
4( J + M) < A' < 8( J + M) . 



0, (? = 1, or (by 

(25) 



(26) 



To represent the resuhs it is convenient to set the mag- 
nitude of J equal to unity, or, here, J = — 1. In this case 
M is restricted by 



1 + A78 < M < 1 + Kl\ , 
in which case this value of c gives 

K 



(27) 



(0,0,0) 




AF-II 

(1/2, 1/2, 1/2) 



M 



1 



FIG. 3: Minimum free energy configurations for J = — 1 as a 
function of M and K. The wavevectors q are indicated along 
with the labeling of YN for the phases. 



A(q) = -lOAf-4- 



K^ 



4(Af- 1) 



(28) 



When we compare this result with the value of A(q) which 
we get for q^ — and for q^ = 7r/2, we get the phase 
diagram shown in Fig. El 



B. J>0 

For positive J the minimization is more complicated 
because the second and third terms in Eq. H21(l now com- 
pete. For A(q) to be extremal, its gradient with respect 
to q must vanish, so that 







- s. 



= s. 



-8Mc^ - AJ(cy + c^) + {cy -f c^)K/R 



-8Mcy - U{c^ + c^) + {c^ + Cz)K/R 



-8Mc^ - U{c^ + Cy) + (c:r + Cy)K/R 



(29) 



There are obviously many subcases for the extrema 
and we will not consider equivalent solutions which cor- 
respond to changing the signs of all the Cq 's or permuting 
their subscripts. Thus we have four cases: 



S'x = Sy ^ Sz = Q , Case I 



Sx = Sy = Q , .s^ 7^ , Case II 



(30) 



(31) 



s^ = , Sy^O , s^ ^ , Case III , (32) 



and 



s^ ^ , Sy^Q , s^ 7^ , Case IV . (33) 



1. Case I 

When Sx = Sy = Sz = Q, then c^, Cj,, and Cz can each 
assume the values -1-1 and —1. So we have 



1 , Case la , 
= 1 , Case lb , 



so that 



A = 6Af 4- 12 J - 4^- , q= (0,0,0), 
A = 6Af - 4 J , q = (0, 0, 27r/a) . 



(34) 



(35) 



Case la, q = (0, 0, 0), is the F phase of YN and Case lb, 
q = (27r/a, 0, 0), the AF-I phase of YN. 



2. Case II 

In this case, qx and qy can independently assume the 
values or 27r/a, so that Cx and Cy independently assume 
the values +1 and —1. Then we have 



Cx — Cy = 1 , Case Ila 



and 



Cx = —Cy — 1 , Case lib . 
In Case Ila we have 

A(q) = 2Af(l + 2c2)+4J(l + 2c^c,) 



-2i^V2 + 2czCx , 
so that minimization with respect to qz yields 

2Kcx \ 



= -8A/c^ - 8Jc.r 



(36) 

(37) 

(38) 
. (39) 



V2 + 2czCx j 

So either Sz = (which repeats Case I), or (since c^ = 1) 

2K 

SMczCx + 8J = -^==== . (40) 



\/2 + 2czCx 



This gives 

32M ^c^ + 6UMc:oCz + 32 J^ 



K^ 



(41) 



l + CzCx 

For K = 0, this gives CxCz = —J/M and 

A(q) = 2M + 4J-4jVAf , (42) 

which is the H< 100 > phase of YN with wavevector 
(q^^jO, 0). For K — this has no range of stabihty. For 
iiT ^ we evaluate A(q) by solving Eq. (|41|l numericaUy. 
In Case lib we have 



A(q) = 2Af (1 + 2c^) - 4J 



(43) 



For positive M we thus have q — (0, 2, l)7r/a, which is 
the AF-III phase of YN and 



A(q) = 2Af - 4 J , Case lib 



(44) 



We discard the case when M is negative because it re- 
peats case lb. 



3. Case III 

Here Cx = ±1 (we need only consider Cx — +1) and Cy 
and Cz are nonzero, determined by 

= -8Mcy - 4J{cx + cz) + {cx + Cz)K/R , 

= -Mdcz - 4:J{cx + Cy) + (cx + Cy)K/R . (45) 

Subtracting and adding one equation from the other we 
get 

?>M{Cy-Cz) = Aj{Cy-Cz)-K{Cy-Cz)lR 

8M{cy + Cz) = {2cx+Cy + Cz){-U + 2K/R) ,{46) 



so that 



Cx -L 7 Cy Cz ^ 



Case Ilia 

8Micy + Cz) = i2 + Cy+Cz)i-4J + K/R) , (47) 

Cx = 1 , 8Af = 4 J - K/R , Case Illb 

8M{cy + Cz) = (2 + Cy + Cz){~AJ + K/R) . (48) 

In Case Ilia we have c^^ = 1, Cj, = c^, where 

8Mcy = -4 J(l + Cy)+K , (49) 

so that 

c,^Cy:^{K~ 4J)/(8Af + 4J) , Case Ilia . (50) 

For this case to apply, we must satisfy the restriction 

\K-4:J\ < |8M-|-4J| . (51) 

Then we obtain 

{K - 4 J)2 



A(q) = -2M-2K 



8M + 4 J 



(52) 



For ivT = this solution is II< 110 > of YN. 
In Case Illb we have Cx — I, 

8M = 4J - K/[{1 + cy){l + Cz)]^/^ , (53) 

and Eq. H48|l becomes 

8M{cy + Cz) = -8M{2 + Cy+Cz) . (54) 

This gives M = or 



Cy + Cz = -1 



(55) 



Since Cx + Cy + Cz — 0, which will appear as Case IVc, 
below, we do not consider it further here. 



4. Case IV 



Here 



= -8MCx - 4:j{Cy + Cz) + (Cy + Cz)K/R , 

= -8Mcy-AJ{cx + Cz) + icx + Cz)K/R, 

= -8Mcz - U(cx + Cy) + {cx + Cy)K/R . (56) 



This set of equations is of the form 




= 0, 



(57) 



where A = -8M and B = -4J + K/R. Note that the 
eigenvalues of this matrix are A + 2B, A — B, and A— B. 
The solution of this set of equations is either of type a (in 
which Cx = Cy = Cz = 0) , type b (in which the eigenvalue 
A + 2B is zero), or type c, (in which the eigenvalue A — B 
is zero). 

The solution of type a is case IVa, with 



Cx = Cy = Cz = , A(q) = -6Af - 2K 



(58) 



For K ^0, this is AF-11 of YN. 

The solution to Eq. H57|l of type b is Case IVb with 



{cx,Cy,Cz) = {c,c,c) , ^ + 2B = 
Setting A + 2B = leads to 

-,1/2 



K 



4( J + M) 



1 



1/2 



and wc have the constraint 

K 



1< 



4( J + M) 



<2, 



(59) 



(60) 



(61) 



and therefore this regime does not appear in the limit 
studied by YN. Then 



A(q) 



-10Af-4J 



K' 



4( J + M) 



(62) 



The solution to Eq. H57() of type c requires A — B = 
and the solution must be a linear combination of the 
two associated eigenvectors. So we introduce Potts-like 
variables^ 



\^X^ Cy , Cz ) 



Note that we have 



/? 



(3 2/3 



V2 V6 ' %/2 ^/6 ' %/6 



2 I 2 I 2 
^x ''' ^y ' ^z 



a'+P' 



and 



CxCy + (cx + Cy)cz = --(a^ + /3^) 



The equation A — B = is 
= 4 J - 8A/ - 

If we write 



K 



v/l-(a2+/32)/2 



4 J = 8M + ^X , 
where ^ can not be negative, then 

a^ + /32 = 2 - 2^-2 = X^ 
This indicates that ^ > 1 or 

K 



0< 



< 1 



Thus we set 



4 J - 8M 



X cos 6* , /3 = X sin ( 



.(63) 

(64) 
(65) 

(66) 

(67) 
(68) 

(69) 

(70) 



where the restriction on 9 will be discussed. These eval- 
uations give 



X^ = 2 - 2 


K 


z 
= 2 - 2{K/K,f , (71) 


so that 


\(ri\ OM 4 7 






J~8M' ^^^' ^^' ^^^' 



and wc have the constraint of Eq. (|69|) . Then, Eq. I|63|) 
gives 

^ ^ \^x 1 ^y^ ^z } 

2X 



V6 



sin(6'-27r/3),sin(6l + 27r/3),sin6'] . (73) 



Note that 9 is arbitrary, so this minimum is realized 
along a curve in wavevector space. Equation H68|l shows 
that X < \/2. As long as X is less than •\/6/2 (i. e. 
Kc > K > Kc/2), all values of 9 are acceptable. If X lies 
between ■\/6/2 and V2 {i- e. Q < K < Kc), then values 
of 9 symmetric around 9 = (and also around 9 = A:7r/3, 
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FIG. 4: Allowed and forbidden (indicated by an F) regions 
9 (in degrees) versus K/ Kc, where Kc = 4 J — 8A/. 



where k is an integer) in which none of the c's exceed one 
in magnitude are allowed. The allowed region, —9c < 9 < 
9c, is determined by the condition 



(2/%/6)Xsin(6ic + 7r/3) = 1 



(74) 



The allowed regions of 9 are illustrated in the Fig. 0] We 
will call this phase the M' phase. 

Although 9 is arbitrary, to get the wavevector seen by 
experiment, we want to have 



qa = (tt — 2ttS , vr + 2ttS, tt) 



(75) 



so that 



^x 

Cy 

Cz 



cos(7r/2 — ttS) = sin(7r(5) , 
cos(7r/2 + ttS) = — sin(7r(5) 
cos(7r/2) = . 



This corresponds to 6* = tt in Eq. (|73|l . so that 



sm{7T6) = (2X/\/6) sin(7r/3) = {X/V2) 



1- 



K 



4 J - 8M 



-,1/2 



1- 



K 
'k'c 



(76) 



1/2 

(77) 



Presumably fluctuation effects or further neighbor inter- 
actions select 9 = TT from the degenerate manifold of all 
9 which nrinimize A(q) and we will consider the second 
mechanism in the Sec. III. 



5. Comparison of Extrema 

To find the global minimum of the eigenvalue, we must 
compare the values of the functions at the above extrema. 
For this purpose we summarize the results in Tabled 

Since Case Ila is hard to analyze analytically, we had 
recourse to a computer program to compare the various 



TABLE I: Free energies of the various states. 



Case 


Energy 


q 


For 


YN Case" 


la 


6M + 12 J - 4:K 


(0,0,0) 


Any 


F 


lb 


6M-4J 


(1,0,0) 


Any 


AF-I 


Ila 


Complex'' 


(g:.,o,o) 


\J/M\ < 1" 


H<100> 


lib 


2M-4J 


(0,1,1/2) 


Any 


AF-III 


Ilia 


2M 2K (-^-"jj^ 
zni zj\ SM+4J 


{0,<ly,qy) 




(8M+4J) 


< 1 


H<110> 


IVa 


-6M - 2K 


(1/2,1/2,1/2) 


Any 


AF-II 


IVb 


in If -If ^^ 


{qx,qx,qx) 


1 ^ 4(,/+M) ^ , 

2 < K < 1 




lUJU IJ i(j+M) 


IVc 


"^Af 4 T ^^ 


(i-M+5,ir 


U ^ A 1 QI\/r 




^lU U 4j_8M 






^.„ 



a For /C = 0. 

b For A" = 0, the energy is 4 J + 2M - 4J^/M. 

c This restriction is only for K — 0. 

d For ii' = we have a two-parameter multiphase: Cx + Cy + Cz = 0. For K j^ 0, the single parameter 9 of Eq. I73II is not fixed. 



local extrema and select the global minimum. Having 
done that, we checked some of the results analytically 
with the result shown in Fig. [S] It is interesting that 
turning on K immediately renders the AF-I and AF-III 
phase unstable. Note also that the M' phase (which is 
the one we want for CEAL) does occur for K as large as 
4J and M being quite small. This agrees with the idea 
that the interactions decrease with increasing separation 
so that \K\ > \J\ > M. 

Finally, we remark that the model with only K and J 
nonzero lacks wavevector selection for K < 8J because 
A(q) only depends on R. So this means that c^, Cy, and 
Cz range over a two parameter manifold of fixed R. 




III. FURTHER NEIGHBOR INTERACTIONS 

In this section we consider the effect of third- and 
further-than-fourth- neighbor interactions. 

A. A "4J" MODEL 

We start by including third-neighbor (tnn) interaction 
coefficients L, as shown in Fig. ^ so that we have a 4 J 
model. This interaction occurs at separations equivalent 
to (3, 1, l)a/4. Including them does not affect Jii(q) or 
^22 (q) but now 

Ji2(q) = J2i(q)* =e-"^('^-+'^«+'=-)/'^[ii:$ + i*](78) 

where $ is as before, but now 

* = e(3lT) + e(3Tl) + e(311) + e(3TT) + 7^ , (79) 

where e{lmn) — exp[i(lkx + mky + nkz)a/A\ and V indi- 
cates inclusion of the permutations V = Cx ^ Cy ^ Cz ^ 



FIG. 5: (Color online) Minimum free energy configurations 
for J = 1, as a function of M and K. The components of 
the wavevector (in units of 27r/a) are the triad of numbers in 
parentheses. The point K = Q, M = 1/2 (filled circle, online 
red) is a two-parameter multiphase point where all states sat- 
isfying Cx + Cy + Cz =0 have minimal free energy. For K j^ 0, 
the region (online magenta) labeled M' is the multiphase re- 
gion in which the single parameter 6 can be chosen as in Fig. 
|1| The line segment (online blue) M — 0, K < 8 is a, mul- 
tiphase region in which only the value of 7? = A'/4 [see Eq. 
12111 is fixed. The phase transitions are continuous except for 
those at M = 0. As soon as K is nonzero, the AF-I phase 
gives way to the H< 110 > phase having small q. and the 
AF-III phase gives way to the M' phase which has wavevector 
within the allowed range of 6 of Eq. 11731 . 



Then 



A(k) = 4:j{CxCy + CyC^ + 



-6M 



+AM{cl + cl + cl) - 2KR , (80) 



where now 

R" = [l + A] + {L/K)A+{L/KfU, (81) 
with 

A = i[$** + $**] =-6 + 2A + 4B + 4C (82) 



n = 



_vl/*\l/ 

4 



9-7A-8B + 8C + AD + 8E -AF . (83) 



Also 



A 
B 
C 
D 

E 
F 



2 _j_ 2 I 2 
C^ -\- Cy -t Cj- 

\Cx H" Cy -r Cz jCxCyCz 

(c^ + cy + c,)icl + cl + cl) 

I 2 2 . 22| 2 2\ 

(4 + 4 + 4). 



(84) 



We win not pursue the analysis to the same level as 
for the 3J model. Here we will show that for an interior 
point in Cx,Cj,,Cz space (i. e. when all these variables 
are less than 1 in absolute value), there is no extremum 
of A(q) for which all the c's are different from one an- 
other. (Thus this model can not give the state of the 
type c = (c, — c, 0), observed for CEALi^) For this anal- 
ysis we consider the equations dX{q)/dca — under the 
assumption that all the c 's are different from one another. 
For a — X this derivative condition is 







AJ{cy + Cz) +8Mc^- 



R 



K{Cy 



+ L{2Cx + 2Cy + 8Cj: + ACyCz{Cx + Cy + Cz) + AC^CyCz 



+ {L^/K) 



'\&Cx — 7cy — 7cy — 56ci^ -I- IGcxCy + 16cx4 



+ {I2cl + 8CyCz){Cx + Cy+ Cz) + 4c^ -I- 8CxCyCz 



(85) 



We now subtract from this equation that which one gets 
by the permutation V and divide the result by c^: — Cy, a 
quantity which, by assumption, is nonzero. Thereby we 
obtain 



= -4:J + 8M- ^i-K + L[6-4:Cz{' 
R \ 



Cx + Cy + Cz ) 



+ iL^/K)[-9 - 56(4 + CxCy + 4) - WcxCy 

+ 164 + 12(C2; + Cy){Cx +Cy+ Cz) 

-8czicx + Cy + Cz) + 4{cl + 4 + c2) ] . 



(86) 



Now, again subtract from this equation that which one 
gets by the permutation V and divide the result by c^—Cz, 
a quantity which, by assumption, is nonzero. Thereby we 
get 

= L{cx + Cy + Cz)-A8{L''/K){cx + Cy + Cz) (87) 



which indicates that Cx + Cy + Cz 



0. Therefore we 



introduce the Potts representation in the form 
X[sin6l + ^/3cos6'l 



Cii — 



V3 



V6 
X [sin 61 -V3 cos 6* 

C 2Xsin6l 
~ %/3 n/6 



If there is an extremum for which all the c's are differ- 
ent from one another, the calculation we have just done 
shows that it must occur for ^ = 0. Rather than fur- 
ther analyze the derivative conditions, it is instructive to 
consider A(q) in terms of these Potts variables. We have 



A 



e 



1 



B = ^^ - 



D 

E 



e 



-X' , 
2 

2^ 
3^'^X^ - 



3\/2 



£,X^ sin{W) , 



—^iX^ sin(36') 
v2 



42 



^X^ sin(36l) , 



Thus 

A = -6 + 6C^ 
9 - 15^2 - 



n 



2£,^X^ 



-3X^ 
'-X-- 



Ix^ 

2 



2^2 



iX'^ sin(36l) . (89) 



3 

8^ 



C - 2i^X^ 



2^/2 



ix^w\{w) 



2\/2CX3 sin(36l) 



(90) 



The important point is that, correct to leading order in 
L/K we have 



R ^ f{e,X') ~ y^{L/K)^X^sm{3d) 
and therefore a contribution to A(q) of 
2\/2, 



<5A(q) 



-L^X^siniSe) 



(91) 



(92) 



Thus the fact that L is nonzero leads to a nonzero term 
in A(q) which is linear in ^ and renders the manifold 
Cx + Cy + Cz = unstable. Since the quadratic term in ^ 
is of the form 



l<-2 



A(q) ~ 1/2X5-^^ 



(93) 
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with Xe^ ^ 8 J — 2K, we see that, when a minimization 
with respect to ^ is performed, one has 



e = ^^^Lx',i,,(^^e) 



(94) 



and now we generate the following term in A(q) which 
depends on 0: 



A(q) = fix') 



9 



L^X''sin^(36l) . 



(95) 



When L is nonzero, it generates a nonzero value of ^, 
i. e. it would take the extremum slightly out of the 
plane Cx + Cy + Cz = 0. But, according to Eq. (|57|l . the 
minimum with the c^'s being unequal, can only occur 
in the plane Cx + Cy + Cz — 0. So, if a minimum can 
not occur in this plane, it can not occur anywhere in the 
interior of c space. In addition, even if a small displace- 
ment out of this plane were allowed (and it would not be 
totally unacceptable in view of the experimental data if 
^ were small enough), the ^-dependent term in Eq. H95() 
favors 6 = 7r/2, which would give a wavevector of the 
form aq/{2n) = (1/2 - S, 1/2 - S, 1/2 + 26), which the 
experimental data do not permit. Within the 4 J model 
this problem can not overcome because the sign of this 
anisotropy in 9 can not be adjusted (it enters in terms of 
positive definite quantities). 



B. STILL FURTHER NEIGHBOR 
INTERACTIONS 

The preceding calculation, although unsuccessful in 
producing an explanation of the data, is nevertheless in- 
structive. It indicates that we need to focus on the sixth 
order anisotropy (in c space) coming from further neigh- 
bor interactions. This requires a term involving six pow- 
ers of the c's. The leading candidate for such a term is 
the exchange interaction at the separation (a, a, a). From 
these eight equivalent neighbors, with exchange constant 
Q, one finds the additional contribution to J„„(q) to be 

(5Jrm(q) = SQ C08{aqx) C08{aqy) cos{aqz) 

= 8Q(2c2-l)(24-l)(2c2-l), (96) 

which leads to an additional term in A(q) whose depen- 
dence on 9 is of the form 



SX{q) 



MQclclcl 



= Q\ 

32 

27 



8^^ 8^X2 



^3^/3 2V3 3%/6 
QX^sm'{ie) . 



X^sin(36i) 



(97) 



(Here we dropped the less significant terms proportional 
to sin(36').) The sign of this term is adjustable and it 
will have the opposite sign from the anisotropy due to L 



TABLE II: Values of the exchange integrals (for J = 1) which 
give values of c close to the observed values of c. Note that a 
small change in Q causes the anisotropy in c-space to change. 



\K\ L M Q 


Cx Cy Cz 


3.50" -0.04 0.016 0.0010 
3.50 -0.04 0.016 -0.0010 
3.50 -0.04 0.017 0.0010 
3.50 -0.03 0.016 0.0010 


-0.339 0.001 0.338 
-0.425 0.215 0.215 
-0.333 0.000 0.333 
-0.360 0.000 0.360 


3.50=^ -0.08 -0.012 0.0050 


-0.337 0.000 0.337 


3.00 0.08 0.131 0.0010 
3.00 0.08 0.131 0.0006 
3.00 0.08 0.131 0.0000 


-0.344 0.022 0.321 
-0.358 0.033 0.324 
-0.410 0.203 0.206 


3.00" 0.00 0.102 0.0001 


-0.333 0.000 0.333 


3.00" 0.04 0.116 0.0004 


-0.340 0.000 0.340 


3.00" -0.04 0.085 0.0003 


-0.336 0.000 0.336 


3.00" -0.08 0.066 0.0012 


-0.338 0.000 0.338 


3.00" -0.12 0.046 0.0025 


-0.338 0.000 0.338 


2.50 -0.08 0.160 0.0008 
2.50" -0.08 0.134 0.0008 


0.000 0.000 0.000 
-0.335 0.000 0.335 


2.50" -0.12 0.114 0.0020 
2.50 -0.12 0.150 0.0020 


-0.338 0.000 0.338 
0.000 0.000 0.000 


2.00" -0.12 0.180 0.0020 


-0.342 0.001 0.341 


2.00" -0.08 0.200 0.0008 


-0.338 0.000 0.338 


1.50" -0.08 0.268 0.0004 


-0.334 0.000 0.334 


1.50" -0.12 0.247 0.0020 


-0.336 0.001 0.335 



a) For this line of parameters, adding —0.002 to Q takes c 

from the (c, — c, 0) phase into the (2c, — c, — c) phase (e. g. 

see the first and second lines of this table) . 



if Q is positive (antiferromagnetic) . It will then favor 9 = 
n7r/3 in Eq. H88|) and if this anisotropy dominates, then 
the wavevector will be of the desired form: aq/(27r) — 
(1/2 -<5, 1/2 + ,5, 1/2). 

If we only consider L and Q, then the condition that 
this anisotropy have the correct sign to explain the 
wavevectors of CEAL is that 



g > (3/8)x5i' 



(98) 



Since the model now has four parameters, we did not 
pursue a definitive numerical analysis of the minima. 
However, to corroborate this argument is sound, we give 
in Table ^] some values of the input parameters which 
give c = (c, — c, 0), for c equal to the experimental value 
c = siuTrJ = 0.338 for S = 0.11,^ This table illustrates 
the phase transition in the anisotropy in c-space which 
takes place for small Q at a value close to that predicted 
by the approximate bound of Eq. (|98|l . Note that for 
small values of L and Q, the values of the other parame- 
ters which give the desired form of c can not be far from 
the region M' of Fig. |31 
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C. EXPERIMENTAL DETERMINATION OF 
PARAMETERS 

In principle we can fix the magnitudes of tfie dominant 
exchange integrals by relating them to several experimen- 
tally observed quantities. These quantities include the 
value of the ordering temperature, Tc^ the Curie- Weiss 
temperature, 0c-Wi for the susceptibility 



X -- C/(T - Gc-w) 



(99) 



and the high-temperature (compared to Tc) specific 
heatiSii22i We consider these in turn and will obtain an es- 
timate for the largest exchange constants J and K (which 
here we denote Kf) to avoid confusion with the symbol 
for Kelvin temperature units.) Crudely we estimate that 
due to fluctuations not included within mean field theory 
the actual ordering temperature, 3.8K, is about 2rMF/3, 
so that Tmf ~ 6K. From Eq. ({T^ we deduce that (ne- 
glecting L and Q) 



24K = 4TMF = -A(q)/fc 
= (-1.88i^o + 0.46J 



5.10M)/fc , (100) 



where we took Kq to be negative (bearing in mind the 
discussion of Fig. [SJ, and we evaluated the constants for 
c = (0.338,-0.338,0). If we only take into account the 
interaction Kq, we get K^/k — —13 K. To see what zero- 
temperature splitting, A_E, of the doublet this implies, 
note that both Tmf and Ai? are proportional to J(q), 
the Fourier transform of the exchange integral. This type 
of relation leads to 



A^ 



S+l 



-m 



MF 



2kT^ 



MF 



(101) 



so that lS.E/k — 12K. This nearly agrees with the result 
AE/k « 15K, given by Boucherle and Schweizeri^i 

Next we consider the Curie- Weiss temperature. This 
is a particularly good quantity to compare to calcula- 
tions because, being the first nontrivial term in the high 
temperature expansion of the uniform susceptibility, it 
is not subject to fluctuation corrections. In Appendix 
B we give a generalization of Eq. H99|) which takes the 
crystal field splitting into account. There we show that 
the Curie- Weiss intercept extrapolated from values of the 
susceptibility x at infinite temperature is related to the 
exchange constants via 



^J„/fc = (20/21)9 



c-w 



(102) 



Following reference ISa we set the Curie- Weiss intercept 
equal to -33K. But, as shown in Appendix B, to get this 
value when an is made from data at T < 300K (rather 
than from infinite temperature), it is necessary to take 

- 28.5K = Y^ J^j/k w [AKo + 12J)/fc . (103) 



If we neglect M, then Eqs. 
determination 



((TUni) and itTU^ lead to the 



Ko/k = -11. 3K , J/k = 6. IK 



(104) 



The value of Kq is fixed to within about 10% by Eq. 
(|100|l . but the value of J is subject to larger (say 20%) un- 
certainty. A question which we can not settle is whether 
it is justified to rely on a pure Heisenberg model to inter- 
pret that Curie- Weiss susceptibility. Attributing contri- 
butions to the susceptibility to the conduction electrons 
or to the diamagnetism of core electrons would somewhat 
modify our estimates. 

The magnetic specific heat C for a system governed by 
the spin Hamiltonian TL gives rise to the limiting value 
CT^ /k at infinite temperature given by 



CT^/k = 



TtH^ l^Trj2(S, .S,)2 



Trl 



Trl 



^E4[^(^+i)]^ = f E^ 



2 

32 ^^'^ 

I, J j 

3N 

= -^[4:K^ + 12 J^ + 6M2 + UL^] , (105) 

where N is the total number of Ce ions. This quantity 
might not be easy to determine experimentally because 
it requires separating off from the total measured specific 
heat (in the temperature range, say, lOK < T < 20K), 
the amount attributed to the lattice and conduction elec- 
trons. 

Finally, we should mention that the interactions we de- 
termine are those renormalized by virtual excitation to 
excited crystal field states. Normally, one might ignore 
such effects. However, as we show in Appendix C, the 
contribution to J from these virtual processes is of the 
same order as we have just determined by our fit to exper- 
iment. These virtual process also imply that long range 
interactions must be present even if one does not invoke 
RKKY interactions. So our appeal to the Q interaction 
[at separation (a, a, a)] is not unreasonable. 



D. APPLICATION TO TmS 

At this point we recall that wavevector selection in 
TmS is of the same form as for CEAL [see Eq. (^], but 
with S — 0.0753 In TmS the Tm spins form an fee lat- 
tice, so the lattice geometry is not the same as for CEAL 
and for TmS the interactions K and L do not occur. 
However, TmS is similar to CEAL in that one can imag- 
ine the dominant exchange interactions limiting one to 
be close to the subspace Cx + Cy + Cz = 0, in which case 
a major concern is to have the anisotropy in wavevector 
space, as in Eq. (|97|l . so that the incommensuration is 
of the form c — {S, —6, 0) rather than c = {6, S, —2S). We 
illustrate this analogy by a brief numerical survey of the 
selected wavevector as a function of the interaction Q (for 
separation (1, 1, 1)) as in Eq. H9t)|l . The result in Table 
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TABLE III: Wavevectors at which A(q) is minimal for values 
of the exchange interactions (in arbitrary units) for the listed 
separations (in units of the lattice constant) on an fee lattice 
as a function of the (1, 1, 1) interaction analogous to Q of Eq. 
1951 1. All interactions are positive (antiferromagnetic). 



Exchange Interactions 


cos{qaa/2) 


(i,i,0) (1,0,0) (i,i,l) (1,1,0) (1,1,1) 


Cx Cy Cz 


2.000 1.000 0.150 0.044 0.003 


-0.304 0.000 0.304 


2.000 1.000 0.150 0.044 0.002 


-0.334 -0.002 0.336 


2.000 1.000 0.150 0.044 0.001 


-0.412 0.150 0.264 


2.000 1.000 0.150 0.044 0.000 


-0.444 0.224 0.224 



This model is somewhat misleading in that it has much 
higher symmetry than that required by crystal symme- 
try. When more general interactions are present, the 
eigenvector of the quadratic free energy matrix associ- 
ated with the eigenvalue which first becomes nonpositive 
as the temperature is lowered determines the form and 
symmetry of the long range order. This critical eigenvec- 
tor must transform according to an irreducible represen- 
tation (irrep) of the symmetry group of the crystal, as 
is discussed recently by one of usi2& This discussion tac- 
itly assumes the impossibility of accidental degeneracy 
wherein two or more irreps having different symmetry 
could simultaneously condense. Accordingly we expect 
that 



mil shows again the effect of this term on the anisotropy 
in wavevector space which can be invoked to explain the 
pattern of incommensuration similar to that of CEAL. 
In addition, we mention that like CEAL, no higher har- 
monics, especially at wavevector 3q were detectedi^ We 
propose that, as we show in the next section, this could 
be understood if the magnetic structure of TmS were to 
consist of the superposition of exactly two wavevectors, 
as is the case for CEALi^ 



IV. QUARTIC TERMS IN THE LANDAU FREE 
ENERGY 

In this section we analyze the quartic terms in the Lan- 
dau free energy in order to investigate the coupling be- 
tween wavevectors in the star of q. Before starting this 
complicated calculation, we describe briefly the physical 
effects we will address. As the temperature is lowered 
in the ordered phase, the effect of the quartic terms in 
the Landau free energy, which is to favor fixed length 
spins, progressively increases. This phenomenon is par- 
ticularly significant for incommensurate systems. For 
many systems having uniaxial anisotropy, order first oc- 
curs in which the spins are aligned along the easy axis 
with sinusoidally modulated amplitude. ^° In that case, 
when the temperature is sufficiently lowered so that the 
fourth order terms become important, the fixed length 
constraint causes the appearance of transverse spin or- 
der, which implies a phase transition, ^"^ and Ni3V208 is 
a recent example of this phenomenon. ^^ As we shall see, 
in CEAL the fixed length constraint favors the simultane- 
ous appearance of incommensurate structures at the two 
wavevectors which combine properly to minimize fiuctu- 
ations in the spin lengths. To show this analytically is 
algebraically quite complicated, as will become apparent. 
(If we only wished to show that the double q state was 
favored relative to the single q state, as was done in Ref. 
UtI the calculation would be much simpler. However, our 
aim was to show that the double q state was favored over 
all other possibilities.) 

In the preceding section we discussed wavevector se- 
lection within a model of isotropic exchange interactions. 



S;(R) 



L^ 



'5'a(qn)e 



zq^-R, 



-I- C. C. 



(106) 



where S'^ (R) is the spin vector at the rth site in the unit 
cell at R, c. c. indicates the complex conjugate of the 
preceding terms, and the sum is over the 12 wavevectors 
(\i which, together with — q^, comprise the star of q. For 
some purposes it is convenient to divide the q^'s into 
three classes Q^ for [i = a, /3, 7 such that for n = 1,2,3,4 



Qn J qn-l-4 = Qn , qn+8 = Q 



7 
n 1 



(107) 



where the Q's are listed in Tables |IV| |V| and EH 

Near the ordering temperature 5'J(qi) can be written 
as a temperature-dependent complex-valued amplitude 
Xi{T) times the critical eigenvector m^(qi) normalized 



byEc 



= 1. Then 



12 



S^(R) 



|a;„m^(q„)e'" 



R 



-|- c. c. 



(108) 



Thus the a;„'s are the complex- valued order parameters 
of this system. The result of representation theory for 
CEAL, given in Ref. 133, is that for the wavevector 
qi = (1/2 — 5, 1/2 + S, l/2)(27r/a), the critical eigenvec- 
tor, which gives the spin components of the two sites in 
the unit cell for the irrep which experimentsi^i^ have 
shown to be the active one, is of the form 



m '(qi) = (ae 
m^^(qi) = (- 



i(p 



, ae 



-icp 






P) 



(109) 



where the real-valued parameters a, /3, and (j) depend 
on the interactions but can be determined from experi- 
mental data. The next step in this calculation is to use 
crystal symmetry to relate the eigenvectors for the other 
wavevectors in the star of q to that given in Eq. (|109|l . 
This is done in the Appendix and the results are listed 
in Tables liVllVl and EH 



TABLE IV: Wavefunctions for the Q" wavevectors/"^ 



n 


Q^ 


mi(QS) 


m2(Q^) 


1 


1 

0. 


-5 


^ + S 


1 

2 


ae'^ 


ae-'"* 


/3 


-ae-"^ 


-ae'^ 


-/3 


2 


1 




-S 


i + <5 


1 

2 


ae~'^ 


ae'"^ 


-/3 


Qe'"^ 


ae-'^ 


-P 


3 


1 


-S 


-k-s 


1 

2 


_Q,gi(7r<S-,^) 


Qg»(-n"<S + 0) 


-/3e*''* 


-Qe^"^ 


ae-^'* 


-P 


4 


1 


-S 


-k-3 


1 


_Q,e'{-«+« 


Qgi(7ri-0) 


/3e*'^'' 


-Qe-'"^ 


ae'^ 


P 



a) Wavevectors are given in units of 2-k /a. 
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TABLE V: Wavefunctions for the Q'^ wavevectors 



(a) 



n 


Qf 


mi(Qf) 


m^(Q^) 


1 


1 

2 


1-5 


1 + 5 


/3 


id) 


ae '^ 


-/3 


-ae-"* 


-ae"^ 


2 


1 
2 


\-5 


i + 5 


-/3 


ae ^ 


ae"*- 


-/3 


ae"* 


ae-'"*" 


3 


1 
2 


\-5 


-1-5 


-/3e'"* 


-ae'(-*-« 


ae»(-*+0) 


-/3 


-ae"* 


ae-"*' 


4 


1 
2 


k-s 


-1-5 


/3e'"' 


_Qg!('ri + 0) 


^gi(7r5-0) 


P 


-ae-"* 


ae"*" 



a) Wavevectors are given in units of 2-k ja. 



TABLE VI: Wavefunctions for the Q'' wavevectors 



(a) 



n 


q;^ 


mi(Q;^) 


m2(Q;^) 


1 


i + 5 


1 
2 


1 
2 


-& 


ae-"* 


/3 


ae'"^ 


-ae"^ 


-P 


-ae-"* 


2 


i+5 


1 
2 


1 
2 


-5 


ae"* 


-P 


Qe-"^ 


ae-'"^ 


~P 


ae''* 


3 


-1-5 


1 
2 


1 
2 


-5 


^gi(7r« + 0) 


-/3e"^^ 


_Qg»('ri-0) 


ae-"* 


-P 


-ae''* 


[4 


-1-5 


1 
2 


1 


-5 


ae'(-*-« 


/3e'"^ 


„ae'(-*+'?') 


ae"* 


P 


-ae-"*" 



a) Wavevectors are given in units of 2n/a. 
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We now turn to the calculation. Equation H15|l shows 
that the fourth order terms in the Landau free energy are 



F4 = 7V-,i&fcT^[S-(R).S-(R)]^ 



(110) 



R- 



where & is a constant of order unity (henceforth we set 
bkT =1). In terms of the order parameters Xi the free 
energy per unit cell is 



F 






F4 , 



fill) 



where x~^ — 4fcT + A(qi) when the small perturba- 
tions to the isotropic Heisenberg model are ignored. At 
quadratic order, there is complete isotropy within the 
order parameter space of twelve complex variables. Our 
objective is to find the direction in the space of the x's 
which has the lowest free energy. This direction will in- 
dicate whether condensation (when ordering takes place) 
takes place via a single wavevector or via the simulta- 
neous condensation into more than one wavevector. To 
study this anisotropy, we will consider the subspace 



Ei^»i' = 



(112) 



where we take c — 1, for convenience. We write 



F, = N-^^Y. E E ^^(91)^0^)5^(93)5^(94) 

X exp[i(qi + q2 + qs + q4) • R] 

= T. T. T.^'c.{<ll)Sl{q2)Sl{q^)Sl{q,) 



Gr gi,g2,g3,?4 a/3 
X0G,qi+q2+q3+q4 



(113) 



where the delta function conserves wavevector to within 
a reciprocal lattice vector G. 

We will decompose F^ into terms involving different 
sets of the critical wavevectors q„ (and their negatives) 
and will express the results in terms of the order param- 
eters Xn- We write 



F4 



E'^ 



(114) 



The first set of terms which we consider are those which 
involve only one wavevector q (By this kind of statement 
we always mean q and — q.) which we denote iSi, where 



'^i = En'E 



E<(q^)' 



ml{ci,) 



(115) 



Next we consider terms involving exactly two different 
wavevectors q^ and q-,-. These are of two kinds, which we 



denote 52a and 526. In the first of these we automatically 
conserve wavevector by taking pairs of opposite wavevec- 
tors. This term (which occurs for arbitrarily chosen pairs 
of wavevectors) is 



52. = ^Y.\^.?\^^"ll\\ll 



l<3 




E™a(^*)"*a(-qj) 



(116) 



The second kind of term is one in which 2qi — 2qj is equal 
to a nonzero reciprocal lattice vector, G. This term is 

526 - E^^'^f EE'^2q.-2q,,G 

i^j r G^tO 



Z^ 



mlM' 



'■] E"^^(-qj) 



]ml{ci^)ml{-qj) 



(117) 



Wavevector conservation in these terms is only satisfied 
when the two wavevectors involved are q2n-i and g2n and 
it is exactly this pair of wavefunctions that are coupled 
in the observed "double-q" statei^ 

There are no terms involving exactly three distinct 
wavevectors. The terms involving four wavevectors, de- 
noted 54 m J involve the wavevectors 



54,1 


(qi, 


-q3,q5,q8) 


54,2 


(qi, 


-q3,q6,q7) 


54,3 


(q2, 


-q4,q5,q8) 


54,4 


(q2, 


-q4,q6,q7) 


54,5 


(qi, 


-q2,q7,-q8) 


54,6 


(qi, 


-q2,-q7,q8) , 



(118) 

the negatives of these, and the set of wavevectors ob- 
tained by the permutation Q" -^ Q,^ -^ Q'*' — > Q", 
which amounts to g„ —J (7,1+4. Here and below n + 4 
is interpreted as n — 8 when n -t- 4 is greater than 12. 
We used a computer program to check that the terms we 
have enumerated are the only ones which can appear in 
fourth order. We write out the first of these: 



54,1 = BxiXgXsZs 



E 



m-(Q?) • m-(Q§ 



m^(Q?; 



m^(Qf; 



m-(Q?).m-(Q5^) m-(Q^ 



m-(Qf; 



m-(Q?) • m-(Qf ) m-(Q^)* • m-(Q^ 



(119) 
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TABLE VII: 



Wavefunctions for each wavevector. We list 



n 


qn 


m^(q„) 


m^(q„) 


1 
2 
3 

4 


h-5 
1-5 




1 

2 
1 
2 

1 

2 
1 
2 


1 
1 


1 
1 


1 
-1 


-1 
1 


-1 
1 
1 
1 


-1 

-1 

-1 

1 


5 
6 

7 
8 


1 

2 
1 
2 

j_ 

2 
1 
2 


\-5 


i+5 

\+5 


1 
-1 


1 
1 


1 
1 


1 


-1 
1 


-1 
1 
1 
1 


9 
10 
11 
12 


i+5 
i+5 


1 

2 
1 
2 

1 

2 
j_ 
9 


k-5 


1 
1 


1 
-1 


1 
1 


-1 
1 
1 
1 


1 


-1 

1 
-1 

-1 



We will now treat the case applicable to CEAL when 
ni'^(q„) is parallel to the appropriate (1,1,1) directionfi*^ 
in which case the wavefunctions are those given in Ta- 
ble IVIII Note that whenever a m(Q3 ) or m(Q4 ) ap- 
pears in one of these fourth order terms, then a m(Q4 )* 
or m(Q3 )* also appears. This means that in using the 
wavefunctions, we may replace e*'^ by unity. Also note 
that the wavefunction for Q^ changes sign for n = 1 on 
going from r = 1 to r = 2, whereas the wavefunctions 
for n ^ 1 do not change sign. This means that any term 
which contains an odd number of variables Q" with n — 1 
vanishes when the sum over r is performed. Thus, out of 
those terms listed above, only 54^1 and 54^4 (their nega- 
tive and their cyclically permuted partners) survive the 
sum over r. We will also need (for t — 1 and 5 = 0) 



M,. 



6Xl'^a(q»)™a(qj) 

a 



(120) 



which we list in Table Ivnil 
Thus we have the result 



Si + S2a + S: 



26 



E 



Xi\ + 



11 ^ 
9 ^ 

n— 1 



^2n-lL-^2nJ 






■El 

l2 



„2^2| 



Mt 



2nV^2n-l\ 



(121) 



and 



^54 



TABLE VIIL 


Matrix elements M 


j of Eq 


nm. 


i/j 


1 


2 


3 


4 


5 


6 


7 


8 


9 


10 


11 


12 


1 


3 


1 


-1 


1 


3 


1 


-1 


1 


3 


1 


— 1 


1 


2 


1 


3 


1 


-1 




-1 


-3 


-1 




-1 


1 


3 


3 


-1 


1 


3 


1 


-1 


1 


-1 


-3 


-1 


-3 


-1 


1 


4 


1 


-1 


1 


3 




3 


1 


-1 




-1 


-3 


-1 


5 


3 


1 


-1 


1 




1 


-1 


1 




1 


— 1 


1 


6 


1 


-1 


1 


3 




3 


1 


-1 




-1 


-3 


-1 


7 


-1 


-3 


-1 


1 


-1 


1 


3 


1 


-1 


1 


-1 


3 


8 


1 


-1 


-3 


-1 




-1 


1 


3 




3 


1 


-1 


9 


3 


1 


-1 


1 




1 


-1 


1 




1 


-1 


1 


10 


1 


-1 


-3 


-1 




-1 


1 


3 




3 


1 


-1 


11 


-1 


1 


-1 


-3 


-1 


-3 


-1 


1 


-1 


1 


3 


1 


12 


1 


3 


1 


-1 




-1 


3 


-1 




-1 


1 


3 



+ (Af9,iiMi4 + M9iMii,4 + NhiMii^i)xQxl^xiXi 

+ (Mio, 12*^23 + MlO, 2^12,3 + M io,3Mi2,2)xioa;*2X2X3) 

-l-c. c. 



44 



■ 2 

/ ^ X^n+lX/inj^'^X^Yi+bXin+S, 
n=0 



+2;4n+22^4n+42^4n+62^4n+7 



C. C. 



(122) 



where here and below the index An + k is interpreted as 
4ri -I- fc — 12 if it is greater than 12. We minimize S by 
fixing the phases optimally, i. e. so that 



e-'^'^r^n-i 



e*"/V2„ , 



a;2n-i = e "■■''r2„-i , a;2n = e""' 'r2„ , (123) 
where all the r's are real and nonnegative. Then 



S 



12 

E^f 

j=i 



22 Y^ ^2 2 

~^ 2_^ ' 2n-l' 2r, 
n—1 



E 



2 2 

T-T - 
1 ] 



64 
"¥ 



r\{rl + rl) + rl{ri, + rl^) + rlirl + r\^) 



12/21 2\i 22, 22 , 22 , 22 



rirsr^rs + r2rAr&rj + r5rrrgri2 



+r6rsriorii + rgruriTi + riori2r2r3 
This is to be minimized under the constraint 

12 



(M13M58 + M15M38 + ]VhsM53)xixlx5Xs 



E' 

i=l 



= 1 



(124) 



(125) 



(M24M67 + M26M47 + M27Mie)x2xlxeX7 

(M57M9,i2 + M^gMr^u + M5^i2Mrg)x5X;XQXi2 

{MesMioM + MesoMs.ii + M6,iiM8,io)a;6a;82;ioa;ii 



To do this write S — Sa + Sb , where 

6 / \ 2 

Sa = 2 + ^ r^Vi-r^ 
fc=i ^ 



9 



E^ 

k=l 



2fe-l' 2fc 



16 









(126) 



where the prime on the summation means that we omit 
terms for which i = 2fc — 1 and j = 2k, and 



^.^#E 



9 



n=0 



{r4n+ir4n+5 — Rn+3^4n+8) 



+ {r4n+2r4n+7 - rAn+ATAn+f,) 
2 



20^ 
"9"^„ 



2 2,2 2 

^4n+1^4n+5 + ''4n+2'"4n+7 



I 2 2 _L 2 2 

+''4n+3''4n+8 + ''4n+4''4n+6 



(127) 



We win minimize 5a with respect to the r^'s. For the set 
of Ti's that minimize Sa, it will happen that the nonneg- 
ative quantity Sb is zero. This shows that this set of r^'s 
minimizes S. 

To minimize Sa we handle the constraint by introduc- 
ing a Lagrange parameter 2A. Then the equations which 
locate extrema of Sa, namely dSA/drn — 4Ar„ ~ 0, are 
(for n= 1,2,3,4,5,6) 



4r2„-i 



4r2„ 



11 



11 



A ^^ 

- 12 



A 



^E^^-^ 



9 



fc=i 







= . (128) 



If both r2n-i and r2n are nonzero, then by subtracting 
their equations, one obtains 



(16/9)(r2Vi-r2j=0 



(129) 



Thus r2„_x = ^2n ^^"i since r„ is nonnegative, we set 



r2n-l = r2n = X„ 



(130) 



Now consider X„ and Xm- Add equations for 7'2n_i and 
r2n and subtract those for r2rn~i and r2„i- Thereby one 
obtains 



2 2 

o l^2n-l + ^2n) + o V 2m-l + ''2rrJ 







(131) 



which indicates that X^ — X^. So for all pairs r2n-i, 
r2n both of whose members are nonzero, we may set their 
X's all equal to X , say. In a similar fashion we show that 
for all such pairs which have only one nonzero member 
we may set the nonzero member equal to Y, the same 
for all such singly nonzero pairs. So we characterize the 
minimum as having k pairs of doubly nonzero members, 
each with value X, and / pairs of singly nonzero members 
assuming the value Y. Then we have that 



Sa = 2 + lY'^ 



hx^-,'- 

9 9 



2k{k-l)X'' 



+2klX'Y^ + {l/2){l~l)Y'' 



(132) 



with the constraint 




2kX^ + IY^ = 1 . 


This leads to the result that 


5. = 2 + l(l-2«»)%l«. + ? 


2k{k- 


+2kX\l 2kX^)+^^^^{l 2kX^f 



(133) 



= AX"^ + BX"^ + C , 



(134) 



where 



4fc2 4, 16,,, ^, 32,2 16,,,, ,, 

A = h-fcH fc(fc-l) fc^H P{l~l) 

9 9 ^ ^ 9 9 ^ ^ 



20fcV(90 - 4fc/3 



(135) 



4fc 16/fc 16fc(/- 1) , ,, ,x , ^ 

B = -- + ^ \^ = -20fc/(9/) ,(136) 



and 



If 



---}-^^^- <-) 



- B/{2A) < Xl^, = l/{2k) , (138) 

then the quadratic form is minimized by setting X"^ = 
—B/{2A). Otherwise, the minimum is realized for X = 
Xmax (for which Z = 0). We see that we never have the 
case of Eq. (|138|l because 



20A:/(90 



B 



A 20k^/{9l) - 4fc/3 k-A - k 



1 1 

> 



(139) 



Therefore the minimum occurs for X = j^Tmax and I — 0, 
where 



Sa - 



1 


9lA+18lkB + 36lk^C 
20fc2 - 12fc; - AOk^ + 2 


36/fc2 

1 
36;fc2 


22 1 


¥~ : 


]k ■ 



■ 88fc^; 



(140) 



So we conclude that the minima occur for fc = 1, and 
for only r2n-i and r2„ nonzero, one sees that Sb — 0, 
so that the minima of Sa are indeed the minima of S. 
These minima correspond to exactly what we want: a 
single pair of equal amplitude order parameters of the 
type we hoped for. 

It should also be noted that the phase difference^ be- 
tween the two condensed waves, given by X2n/x2n-i = 
gi7r/2 g^jgQ agrees with the conclusions of Forgan et al^ 
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that the structures of the two incommensurate wavevec- 
tors add in quadrature. In addition, our calculation sup- 
ports their argument that the variation of the magni- 
tude of the spin over the incommensurate wave should 
be minimial. Our calculation also explains why the fixed 
length constraint does not require substantial values of 
higher harmonics, such as S(3q). However, this picture 
can not be totally correct, because the double-q struc- 
ture does not completely eliminate the variation of the 
magnitude of the spin. The spin structure consists of two 
helices of opposite chirality and the ellipticity of these he- 
lices decreases with decreasing temperature, but the ec- 
centricity of the polarization ellipse extrapolated to zero 
temperature^ is too large to be explained by anisotropy 
alone. Probably some, or all, of this eccentricity should 
be explained by Kondo-like behaviori^ 



V. CONCLUSION 

We may summarize our conclusions as follows. 

• For the fee antiferromagnet with first and second 
neighbor interactions we located a previously over- 
looked multiphase point [see Eq. if^ ] at which 
wavevector selection is infinitely degenerate. 

• We have extended the analysis of Yamamoto and 
Nagamiyai^ to determine the minimum free energy 
of magnetic structures of CeAl2 (which is a two 
sublattice fee incommensurate magnet) for a model 
consisting of three shells of isotropic exchange in- 
teractions. The phase diagram in terms of these 
interactions (see Figs. |21and|31) has an incommensu- 
rate phase with a wavevector in a degenerate man- 
ifold which includes the observed incommensurate 
wavevector for CeAl2. 

• We analyzed the effect of third nearest neighbors on 
the degenerate manifold of the three shell model 
and found that it gave the wrong anisotropy in 
wavevector space to explain the data for CeAl2. 
However, the correct sign of the anisotropy (which 
would give a wavevector of the form (1/2 — (5, 1/2 + 
6, 1/2) in units of 27r/a), can be obtained if the in- 
teraction Q of neighbors at separation (a, a, a) ex- 
ceeds a rather small threshold value. Since CeAl2 
is a metal subject to RKKYiS interactions, we sug- 
gest that such an interaction is not unreasonable. 
By way of illustration we give (see Table ^ . some 
explicit values of exchange parameters that will 
give the correct incommensurate wavevectors. 

• By analyzing the form of the fourth order terms 
in the Landau expansion, we show that for the 
wavevectors appropriate to the ordered phase of 
CeAl2, the observed "double-q" stated is favored 
over any other combination of wavevector(s) in the 
star of q. This result is not a common one for a cu- 
bic system. In addition our analysis reproduces the 



relative phase observed^ between the two coupled 
wavevectors. 

• By relating the exchange constants to the Curie- 
Weiss intercept temperature 9c-w of the inverse 
susceptibility and to the ordering temperature, we 
developed the estimates for the nearest neighbor 
ferromagnetic interaction, K/k — —11 ± IK, and 
for the next-nearest neighbor antiferromagnetic in- 
teraction, J/k = 6 ± IK. 

• We also showed (see appendix C) that the exchange 
interactions are significantly renormalized by vir- 
tual crystal field excitations. This effect leads to 
rather long-range exchange interactions. 

• It is possible that our analysis of wavevector se- 
lection can explain the similar incommensurate 
wavevector observed22i24 for the Kondo-like sys- 
tem TmS. Although the anisotropy axis is different 
for TmS than for CEAL, one may speculate that 
the fourth order terms in TmS may give rise to a 
double-q state, although such a state has not yet 
been observed in TmS. 



APPENDIX A: SPIN FUNCTIONS FOR THE 
STAR OF Q 

In this appendix we determine the spin functions for 
the different wavevectors in the star of q, given that for 

aq/(27r) = 1/2 - (5, 1/2 + (5, 1/2 (Al) 

the spin functions for the two sites in the unit cell arf^^^s 



mi(q) ^ [ae"l',ae~"'' , /3] 



mi(-q)' 



m2(q) = [-ae-*^-ae^* ,-f3] ,-m2(-q)* (A2) 

where a, /3, and ip are real valued constants, are fixed 
by the interactions through the quadratic terms in the 
free energy. Since we will study the quartic terms which 
couple different wavevectors, we need to tabulate the spin 
functions for the different wavevectors. 

The star of the wavevector consists of 24 vectors which 
are ±Q^, ±Q^, and ±Q;^, for n = 1,2,3,4. These Q's 
are listed in Tables IIVI Ivl and IVII The spin functions 
for different wavevectors are related by the symmetry 
operations of the crystal, which is space group #227, 
Fd3m, in the International Tables for Crystallography 
(ITC).25 

In Eq. IJA2|I we gave the spin wavefunction for Q". We 
now consider the effect on this function of the operation 
(x, y, z) -^ {—y, —X, z) {4l^i7 in ITC), which we regard as 
a mirror which interchanges x and y followed by a two- 
fold rotation about z. Because spin is a pseudovector this 
operation on spin is 



{mx,my,m^) -^ {■my,mx,-mz) 



(A3) 
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Thus, before transformation we have 

mj;(Ri,ri) = 2acos(qj-Rj4 

my(Ri,Ti) = 2acos(qi-Ri- 

inz(Ri,Ti) = 2/?cos(q,; • R,;) 

ma:(Ri,T2) = -2QC0s(qi-Ri 

my(Ri,T2) = -2acos(qi-Ri 

inz(Ri,T2) = -2/3cos(qi • Ri 



(A4) 



where R^ = (X^, Yi, Zi) specifics the location of the unit 
cell before transformation, q; is the wavevector before 
transformation, given in Eq. IJA1|I . and 



XI = (0,0,0), X2 = a(l,l,l)/4 



(A5) 



After transformation (indicated by primes) Eq. (|A3|I 
gives 

m^(R/,rij) = 2acos(qi • Ri — (/)) 

my{Rf,Tij) = 2a cos(qi • Ri + (/)) 

m'^(R/,rij) = -2/3cos(qi • Ri) 

m^(R/,T2,/) = -2acos(qi • Ri + (/)) 

toJ^(R/,T2,/) = -2acos(qi • Ri - (/)) 

m;(R/,T2j) = 2/3cos(q, -R,) . (A6) 

where R/ = (-^/, ^f, Zf). If the initial position is 

r = {x„ y„ z,) + n = (x„ r„ z,) , (A7) 

then the final position is 

r' - (-r„ -X.„ Z,) EE {Xf,Yf, Zf) + Tf , (A8) 

so that Tif = Ti. We now express q^ • R^ in terms of the 
final coordinates: 

qixXi + qiyYi + q^^^Zi = q^^Yf ~ q^yXj + qi^Zj (,A9) 

which can be written as q, • Rj = q/ • R/ , where we have 
(to within a reciprocal lattice vector) 



q/ 



( 1i,y: qi,x,qi,z) — Q2 



(AlO) 



Thus 



m^(R/,ri) = 2acos(Q^R/-0) 
m^(R/,ri) = 2a cos(Q^ • R/ + </.) 
m;(R/,ri) = -2/3cos(Q2 -R/) . (All) 

Now consider Ti = T2. Then if the initial position is 

r = {X, + a/4, r, + a/4, Z, + a/4) , (A12) 

then the final position is 



i-Y, - a/4, -X, - a/4, Z, + a/4) 
iXf,Yf,Zf)+Tf, 



so that, in this case, Tf = T2. 

X/ = -Y, - a/2 , Yf = -X, - a/2 , Zf^Z, (A14) 

We express q^ • R^ in terms of the final coordinates: 

q,; • R,; = qix{-Yf - a/2) + qiy{-Xf - a/2)+qi:,Zi 
= Q^ • R/ - TT . (A15) 

Thus 

m^(R/,r2) = 2a cos(Q J • R/ + (/>) 
TOy(R/,r2) = 2acos(Q^-R/-0) 
TO;(R/,r2) = -2/?cos(Q^-R/) . (A16) 

Thus for wavevector Q2 the Fourier component vector 
(which we put into Table llV)l is 

{ae-"^ , ae"^ , -/? ; ae*"^ , ae"*"^ , -[3) . (A17) 

Next we study the effect of the transformation (x, y, z) —^ 
{x + 1/4, —y + 1/4, z + 1/4). Before transformation the 
Fourier coefficients are those of Eq. ljA4p . Since this 
transformation is a mirror operation we have, after trans- 
formation that 

m^(R/, Ti/) = -2acos(qi • R; + (j>) 

my{Rf,Tif) = 2acos(q., • R., - (j)) 

m'^(R/,ri/) = -2/?cos(qj • Ri) 

m^(R/,r2/) = 2a cos(qj • Rj - (/)) 

m^(R/, r2/) = -2acos(qi • Ri + 0) 

m;(R/,r2/) = 2/3cos(q, -R,) . (A18) 

For Ti — 1 the initial position is 

r = R,+Ti = (X„y„ZO, (A19) 

and, using the transformation, the final position is 

r' = (X, + a/4, -K, + a/4, Z, + a/4) 

^ {Xf,Yf,Zf)+nj . (A20) 



Thus Tif = T2 and 

Qi ■ Ri = Q3 • R./ 

Then 



(A21) 



(A13) 



to;(R/,T2) = -2acos(Q;^ •R/ + (/.) 
m;(R/,T2) = 2acos(Q;^-R/-0) 
m^(R/,T2) = -2/3cos(Q^-R/) . (A22) 

Using the transformation on r = {Xi + a/4, Yi + a/4, Zi + 
a/4), we write 

r' == (X, + a/2, -r„ Z, + a/2) 

= {X, + a/2, -y,, Z, + a/2) + n , (A23) 



19 



so that T2f = Ti and 

Rf = iX, + a/2,-Y,,Z, + a/2) 
Thus 



so that 



(A24) 



Qj • Rj = q/ • R/ - q/ • (a/2, 0, a/2) 

= Q^ • R/ - TT + 7r(5 , (A25) 



so that 



TO^(R/,ri) = -2acos(Qg^ -Rj-^ + TT^) 
m^(R/,ri) = 2acos(Q;^-R/ + (/i + 7r(5) 
m',(Rf,Ti) = -2/?cos(Q^-R/+7r(5) . (A26) 

Thus for wavevector Q3 the Fourier component vector 
(which we put into Table llV)l is 



-ae"*" , ae""*' , -/3) . (A27) 

Next we study the effect of the transformation 
{x, y, z) -^ {-y + 1/4, X +1/4, z + 1/4) (#16 in ITC). 
Since this transformation is a four- fold screw axis, we 
have, after transformation that 

m^(R/,ri/) — — 2acos(qi • Ri — (/)) 

my(R/,ri/) = 2a cos(qi • Rj + (^) 

m^(R/,ri/) = 2/3cos(qi • Ri) 

m^(R/, T21) = 2a cos(qj • R^ + (/)) 

mJ^(R/,r2/) = -2acos(qi -Ri - 0) 

m',(R/,r2/) = -2^cos(q, -R,) . (A28) 

For Ti = T\, and if R^ = (X^, Yi, Zi), we have 

r' - i-Y, + a/4, X, + a/4, Z, + a/4) , (A29) 

so that Tif — T2 and to within a reciprocal lattice vector 
this gives 



q/ = Q 



4 ; 



(A30) 



so that 



TO^(R/,r2) = -2acos(Q^ -R/-^) 
m;^(R/,r2) = 2acos(Q^-R/ + 0) 
m;(R/,r2) ^ 2/?cos(Q^-R/) . (A31) 

For Ti = T2, r ^ {Xi + a/4, Yi + a/4, Zi + a/4) and 

r' = (-y„ X, + a/2, Z, + a/2) , (A32) 

so that T2f — Ti and 

q^-R., = q/l-Hf -{a/2){q^i + q^i) 

= Q^-R/+7r(-l + <5) . (A33) 



m'A^f^Ti) 
m'y{'B.f,Ti) 
m'zi'^f^Ti) 



-2a cos(Q4 • R/ + + nS) 
2acos(Q2 -Rz-^ + ttJ) 
2/3cos(Q4 •R/ + 7r(5) . (A34) 



Thus we have the results for the order parameter wave- 
functions given in Table HVl To get the wavefunctions for 
Q^ and for Q^ is much easier: one simply uses the three- 
fold rotation axis about (111) to get the results given in 
Tables El and ED 



APPENDIX B: CURIE- WEISS SUSCEPTIBILITY 
IN A CRYSTAL FIELD 

Here we develop a formula for the susceptibility correct 
to leading order in the exchange interactions, Jy. For 
this purpose we write the Hamiltonian as 



n 



Ti-o + Xy ^ JijSoS,i ■ Soff, j 






(Bl) 



where the ScS.i is the effective spin 1/2 operator we have 
used throughout our calculations and A, a scale factor 
for the perturbation, is set equal to unity in the final re- 
sults. Here TCq includes all terms for Jij — 0. Thus Ho 
is the Hamiltonian for spins subject to the cubic crystal 
field and the external magnetic field, but with no ex- 
change interactions between neighboring spins. It will 
be convenient to express this Hamiltonian in terms of 
the magnetic moment operator fj,i for site i. We write 
Scff = (3/5)J — [3/(5.gj/i_B)]/x, so that (with gj = 6/7) 

n = no + {7/lo^lBfxY,Jv^^^ ■ f^j , (B2) 

Correct to leading order in A we use thermodynamic per- 
turbation theorji^ to write the free energy as 



F(A) = F(A = 0) + 



49A 
100^1 



^J.,(M,)o(/Xj-)o,(B3) 



i<3 



where F{X = 0) is the free energy for the Hamiltonian 
Ho and 

{X)o = Tr[Ae-'5^"]/Tr[e-^^»] . (B4) 

Then the susceptibility per spin, x = d{^i ) /dH\jj^f^, is 



X(A) 



N- 



dH^ -X(A = 0) 
[49A/(100Ai|)] Y, J,,x(A = 0)2 , (B5) 



where N is the total number of Ce ions. Thus 



x(A)-i = x(A = 0)-i + [49A/(100Aii)] ^^ J,, 



0{\') 



X(A = 0)-i + [49/(100/4)] (4X + 12J + . . .) 
+0{Jl) . (B6) 



To obtain x(A = 0) we took the wavefunctions of the 
ground doublet in the cubic crystal field to be 



200 



|0) = 



/576|5/2, ±3/2) - v^|5/2, t5/2) , (B7) 



in the J, J^ representation. The remaining states form 
the four-fold degenerate excited state at a relative energy 
which we denote kTn. Then we found that 



kx{X = 0) 



Ai| 



320(1 - e-^o/^) 

50 + 260e-^«/^ 

49T[2 + 4e-'^«/'^] 



+ 



(B8) 



At low temperature, the second term displays the Curie- 
like 1/T dependence corresponding to the moment in the 
ground doublet and the first term is the so-called Van 
Vleck temperature-independent susceptibility,"*^ To illus- 
trate the effect of this term, we show the inverse suscepti- 
bility in Fig. Elfor the Ce ion (J = 5/2) in a cubic crystal 
field with a doublet-quartet energy splitting of kTg with 
Tq = lOOK. In the high temperature limit T » Tq we 
have 



X(A = 0) 



9p{J+l)^il _ IB^ij 



3kT 

in which case for A = 1 we have 
7k 



7kT 



(B9) 



with 



15/i| 



ec- 



c-w 



[T - Gc-w] + O { Jy /T, {Tg/Tf} GBIO) 



^21 Y^ Jy /20 w -21{K + 3J)/5 .(Bll) 



In Fig. El we also show the inverse susceptibility when 
8c-w = — 29.9K, a value which gives the Curie- Weiss 
intercept (extrapolated from T = 300) of — 33K as in 
Ref. IM 



APPENDIX C: EFFECTIVE INTERACTIONS 
VIA EXCITED QUARTET STATES 

Here we consider effective interactions which occur via 
excited virtual crystal field states. In a general formu- 
lation one considers manifolds Mn in which n spins are 
in their excited quartet crystal field level, whose energy 
is kTg relative to the crystal field ground state. We are 
interested in the effective Hamiltonian Hq for Mq at low 
temperature and here we discuss its evaluation within 
low order perturbation theory. We define 



Hn 



VnUV^ 



(CI) 



where Vn is the projection operators for the manifold 
M.n- Clearly the lowest approximation is to neglect en- 
tirely all processes except those within the manifold M.q 




100 150 

T(K) 



300 



FIG. 6: (Color online) Inverse susceptibility, 1/x- The dots 
are data points taken from the more extensive data set of Ref. 
35. The lower (online blue) curve is for A = and the upper 
(online red) one is correct to first order in A [using Eq. 1B6H 
for ^ . Jij /k = 28.5K. For this curve the intercept extrap- 
olated from 280K < T < 300 is -33K, as indicated by the 
arrow. The intercept extrapolated from infinite temperature 
is 6'c-w = -(21/20) Y,. Jij = -30K. So even at T = 300K 

there is still a noticeable departure [of order (Tq/T)^] from 
the infinite temperature behavior. 



and this was an implicit assumption of our calculations 
in the body of this paper. Processes involving virtual 
state in the manifold Mi enter via second order per- 
turbation theory. These terms are obtained just as for 
superconductivity^^ with the result that 



^° - ^'''-klhT. 



n Ti-O^nTi-nfi 



(C2) 



One can use this formalism to reproduce the formula for 
the zero-temperature Van Vleck susceptibility tensor;^ 
xjj i. However, our present aim is rather to analyze ef- 
fective exchange interactions which arise in this way. To 
illustrate the phenomenon, consider contributions to Eq. 
(|C2|) when Tio,n is taken to be the nn exchange interac- 
tion. Although we wrote this interaction as XoSoff.i-Soffj 
it really should be represented as K'Ji ■ 3j , where, since 
Soff.i = [{gj - l)/.9o]Ji = 3/5Ji one has that K' = 
(3/5)^1^0 = (9/25)ii'o. Then we obtain a contribution l^y 
to the effective nnn exchange interaction between spins 
i and j [at separation (a/2, a/2,0)] using the nn interac- 
tions Jik between spins i and k and Jkj between k and j . 
Since for a nnn pair i,j there is only one choice for the 
intermediate site k to be an nn of both sites i and j, we 



21 



have 



V^, 



K'' 



kT, 



Q 



■VoJ^ ■ JkVlJk ■ JjVo . (C3) 



Because of the cubic symmetry of the crystal field one 
has 



VoJk,o.Vi.h,pVo = (20/9)(5„,;37'o 



(C4) 



To obtain this result it is convenient to take a = j3 = z 
and use the ground state wavefunctions of Eq. (|B7p . 
Then Eq. ^^ yields 



_ 20(9Xo/25)^ 

-boff.i ■ OeS,\ ■ 



5kTo 



(C5) 



This means that due to these processes we have that 



J^ J 



4Kl 
5kTQ 



= J-5J. 



(C6) 



For Kf:,/k = lOK and Tq = lOOK, this gives 5J/k w 
0.8K. The nn interaction is renormalized in a similar way, 



except that if sites i and j are nn's, then there are now 
two choices for the site k to be an nn of both i and j. 
Thus 



Ko^Ko 



8Kl 
5kT, 



Ko - SKo 



Q 



(C7) 



with SKq = 1.6K. As a final example, we similarly find 
for sites separated by (a, 0, 0) that there are four inter- 
mediate paths of sites separated by (a/2, a/2, 0), so that 



M 



16/2 
M - -—— = M - SM . 



5kT, 



Q 



(C8) 



Taking J/k = 5K, we find that 5M/k « 0.8K, so that 
5M/ J — 0.16, a value which is comparable to those used 
in Table HJ 

These results imply that even if the bare Hamilto- 
nian only has nn interactions, virtual processes involv- 
ing higher crystal field states will induce nnn interac- 
tions approximately of the size we will deduce from fit- 
ting experiments. This mechanism in higher order will 
produce significant longer range interactions even if the 
bare Hamiltonian has only nn interactions initially. 
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